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I.  INTRODUCTION 


^A  program  of  continuing  development  of  several  fluid  dynamic  models  has 
been  carried  out  for  NRL  by  JAYCOR.  Specifically,  major  developments  and 
modifications  of  several  large  computer  codes  have  been  accomplished.  This 
has  been  done  to  provide  NRL  with  an  expanded  capability  of  modelling 
diverse  physical  phenomena,  all  of  which  are  amenable  to  analysis  utilizing 
advanced  hydrodynamic  computational  techniques.^ 

Among  the  program  modifications  performed  in  this  work  are  the  following: 

1.  Inclusion  and  testing  of  thermal  conduction  routines  that  have 
multi -material  capabilities  and  that  utilize  classical  plasma  thermal 
conduction. 

2.  Multi-material  capabilities  have  been  added  to  FAST20. 

/ 

3.  A laser  energy  deposition  algorithm  has  been  developed  and  coupled 
with  both  FAST2D  and  ADINC. 

4.  A new  initializer  that  generates  semi-analytical  density,  temperature, 
and  velocity  profiles  has  been  incorporated  into  both  ADINC  and  FAST2D. 

5.  An  improved  diagnostic  that  allows  for  a much  better  determination 
of  the  Rayleigh  Taylor  instability  growth  rates  has  been  added  to  FAST2D. 

6.  Improved  modelling  of  the  planetary  boundary  layer  has  been  incorpo- 
rated into  both  a 2D  and  NRL's  3D  tropical  cyclone  models. 

Utilizing  these  and  other  expanded  computational  capabilities,  the 
following  studies  have  been  undertaken: 

1.  The  thermal  burn-through  of  multilayered  target  foils  has  been 
studied  for  different  materials  and  thicknesses. 


2.  Rayleigh  Taylor  Instability  analyses  have  been  performed.  Growth 
rates  have  been  determined  and  dispersion  relations  deduced.  The  sensitivity 
of  these  results  to  improved  spatial  resolution  has  been  examined. 

3.  Dynamic  stabilization  of  the  Rayleigh-Taylor  instabilities  arising 
In  laser  fusion  has  been  modelled  and  demonstrated. 

4.  The  instabilities  arising  from  the  use  of  support  webs  in  pellet 
designs  and  those  arising  from  asynmetries  in  the  laser  intensity  have  been 
studied. 


5.  Turbulent  boundary  layer  phenomena  have  been  modelled  and  incorporated 
into  a dynamic  atmospheric  model. 

In  the  following  sections,  both  these  computational  advances  and 
scientific  results  will  be  discussed  in  detail. 
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II.  COOE  DEVELOPMENT 

Both  a ID  and  a 2D  thermal  diffusion  routine  that  utilized  the 
plasma  thermal  conductivity  had  been  developed  at  NRL.  However,  several 
modifications  were  necessary  for  their  more  general  use  In  laser-fusion 
studies.  In  order  to  study  multi -material  fusion  target  designs,  space 
dependence  of  the  thermal  conduction  coefficient,  multiple  density  arrays, 
and  more  convenient  temperature  units  had  to  be  Introduced.  Also,  the 
thermal  conductivity  of  a mixture  of  Ions  had  to  be  calculated  (since  plastic 
and  glass  are  comnonly  used  for  target  construction),  and  the  charge  de- 
pendence had  to  be  Introduced  Into  the  thermal  coefficient  (only  the  trivial 
case  of  D-T  with  Z*1  had  been  considered). 

Another  important  change  needed  was  the  generalization  of  the  boundary 
conditions.  For  several  test  calculations,  It  was  necessary  to  specify  the 
boundary  temperatures.  However,  for  the  laser  fusion  studies,  no  heat  flow 
across  the  boundaries  was  required.  In  order  to  implement  these  different 
boundary  conditions  (and  any  combination  thereof),  the  following  changes  were 
made. 


The  thermal  diffusion  routines  use  a tridiagonal  solver  to  update  the 
temperatures.  Namely,  for  ID, 


Vi  - 1 


+ B1T1 


+ C1T1  + 1 


(1) 


Is  solved,  where  the  T's  are  the  new  temperatures  and  the  coefficients  are 
functions  of  the  old  temperatures  and  other  factors.  At  the  left  boundary, 
for  example,  we  would  have: 

A1T0  * B1T1  ♦ C1T2  • °1  • (2)  -i 

where  TQ  Is  some  to-be-specified  boundary  temperature.  If  we  write  it  as:  Jj 

T0  » a TL  ♦ 6 Ti  , (3) 


then  we  can  set  It  to  some  specified  boundary  temperature  (ot  * 1,  6*0), 
or  assure  that  there  is  no  heat  loss  across  the  boundary  (a  ■ 0,  6 * 1). 


J 
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Substituting  equation  (3)  into  (2)  gives  for  the  new  coefficients: 


B-j  + BA1 


D-j  - a T^A-j 


Of  course  the  same  thing  can  be  done  for  the  right  boundary  with  a different 
a and  6. 

Another  important  modification  that  involved  the  thermal  conduction 
routines  was  the  introduction  of  the  laser  energy  directly  into  the  calcula- 
tions. Previously,  the  laser  deposition  had  been  modelled  by  driving  the 
boundary  temperature  of  the  system.  This  is  inadequate  for  any  detailed 
analyses  of  laser  fusion  experiments  since  all  the  details  of  the  laser 
deposition  are  missing.  By  directly  inputting  the  laser  energy  into  the 
system,  both  the  temporal  and  spatial  distribution  of  the  flux  could  be 
modelled,  and  the  amount  of  energy  deposited  could  be  controlled. 

We  assume  that  all  the  laser  energy  is  absorbed  through  resonant  absorp- 
tion. Hence,  the  critical  density  (where  the  laser  frequency  equals  the 
plasma  frequency)  is  located  and  the  energy  is  deposited  into  the  two  cells 
about  this  critical  point.  This  is  done  by  introducing  an  energy  source 
term  into  the  thermal  diffusion  routines. 

A routine  was  written  that  generated  the  spatial  and  temporal  distribu- 
tion of  the  flux.  This  was  done  quite  generally  so  that  several  different 
temporal  evolutions  and  spatial  illuminations  could  be  employed.  Also,  a 
special  smoothing  algorithm  was  introduced  so  that  when  the  critical 
density  point  crossed  a cell  boundary,  the  finite  resolution  of  the  mesh 
didn't  cause  any  thermal  shocks  to  be  generated.  Finally,  to  simulate 
Inhomogeneities  in  the  laser  illumination,  the  ability  to  introduce  con- 
trolled perturbations  into  the  laser  intensity  was  added. 
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Joseph  Orens  of  NRL  had  developed  a semi-analytical  technique  for 
generating  density,  temperature,  and  velocity  profiles  that  were  quasi- 
static solutions  of  a set  of  hydrodynamical  equations  that  modelled  a 
laser-plasma  interaction.  With  his  help,  subroutines  that  generated  these 
profiles  were  introduced  into  the  initializer  of  both  the  10  and  20  codes. 

In  this  way,  the  codes  could  be  initialized  with  profiles  that  were  known  to 
be  essentially  correct.  The  compatibility  of  these  quasi-static  solutions 
and  the  codes  was  tested.  For  a constant  laser  intensity,  only  minor  changes 
occurred  in  the  initial  profiles  when  they  were  stepped  in  time  in  the  codes 
for  1000  timesteps.  With  these  initial  profiles,  we  could  now  do  our 
Rayleigh-Taylor  studies  with  the  confidence  of  knowing  that  any  large  changes 
in  the  initial  profiles  that  developed  in  time  were  due  to  fluid  instabilities 
which  we  had  initiated  in  a controlled  way. 

Finally,  in  order  to  study  the  growth  of  the  Rayleigh-Taylor  modes  that 
we  initiated,  a better  diagnostic  was  needed  in  FAST2D.  Previously,  either 
the  location  of  the  ablation  layer  or  a time  and  space  averaged  vorticity 
was  used.  For  example,  the  ablation  layer  edge  was  determined  as  a function 
of  y and  then  Fourier  analyzed.  The  Fourier  coefficients  then  determined  the 
mode  amplitudes.  However,  acoustic  waves  were  known  to  be  generated  in  the 
ablation  layer,  and  these  would  shift  the  ablation  edge  enough  to  add  a large 
amount  of  noise  in  the  determination  of  the  mode  amplitudes. 

A much  better  diagnostic  was  found.  The  ablation  edge  was  accurately 
determined  for  each  y and  then  the  mass  contained  up  to  the  ablation  edge  was 
calculated.  This  summed  mass  was  then  Fourier  analyzed  and  found  to  be  much 
less  sensitive  to  the  acoustic  wave  disturbances.  Examples  of  the  quality  of 
this  diagnostic  will  be  given  in  Section  III. 

Two  different  parameterizations  of  the  planetary  boundary  layer  have 
been  incorporated  into  a 2D  tropical  cyclone  model  to  determine  their  rela- 
tive merit.  Since  the  multi-layer  scheme  is  fully  vectorized,  there  is  no 
noticeable  computational  cost  increase  in  its  use.  The  other  parameteriza- 
tion, a matching  algorithm  was  found  to  be  sufficiently  accurate,  and  has 
been  incorporated  into  NRL's  3D  tropical  cyclone  code. 
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III.  SCIENTIFIC  RESULTS 


In  ID,  the  thermal  conductivity  equation  is  written: 

3E  3 v JI 
at  3x  * 3x 


(5) 


where  K,  the  thermal  conductivity,  is  written  as: 
K - KqT5/2  . 


(6) 


For  the  case  of  constant  density  and  zero  initial  temperature,  and  using 
E » 3/2  pT,  Equation  (5)  has  the  solution: 


T » T 

0 xo 


2/5 


(7) 


where  Trt  and  x are  constants  and 
o o 

5'2 

v 'rer 


(8) 


Both  the  ID  and  2D  thermal  diffusion  routines  were  tested  by  driving 
the  boundary  temperature  according  to  Equation  (7)  with  x = 0.  A compari- 
son of  the  analytical  and  numerical  results  for  the  ID  case  are  shown  in 
Figure  1.  The  temperature  profiles  for  several  different  times  are  dis- 
played along  with  the  thermal  "shock"  front  position  as  a function  of  time. 
Note  that  the  agreement  is  excellent. 

Also  the  flux  deposition  algorithm  was  checked  by  depositing  a fixed 
amount  of  energy  Into  the  system  and  then  diffusing  it  with  the  boundary 
condition  of  zero  flux  across  the  boundaries.  Conservation  of  energy  was 
observed  to  one  tenth  of  a percent  accuracy. 

In  order  to  test  the  multi-material  capabilities  of  the  ID  thermal  con- 
duction routine,  the  following  study  was  undertaken.  We  reproduced  the  test 
of  Figure  1 for  a pure  D-T  (deuterium-tritium)  plasma  and  then  repeated  the 
calculations  with  a ly  thick  foil  of  different  materials  Inserted  into  the 


system.  "Shock"  front  position  versus  time  plots  similar  to  the  one  in 
Figure  1 are  displayed  in  Figure  2.  As  can  be  seen,  glass  (S^Oj),  copper, 
and  gold  foils  of  the  densities  stated  were  used.  All  of  these  materials  slow 
down  the  thermal  "shock".  Since  gold  has  the  highest  Z,  ionization,  and 
since  the  thermal  conductivity  is  inversely  proportional  to  Z,  it  inhibits 
the  flow  most  strongly.  However,  at  large  distances  from  the  foils,  all  the 
results  seem  to  display  approximately  the  same  velocity. 

In  order  to  study  the  asymptotic  region  of  Figure  2,  a new  series  of 
runs  was  performed.  The  pure  D-T  results  were  compared  to  runs  where  one  and 
two  micron  slabs  of  glass  were  inserted.  This  comparison  is  presented  in 
Figure  3.  As  expected,  the  thermal  "shock"  front  is  slowed  down  more  when 
passing  through  the  thicker  slab.  At  large  distances  form  the  slab,  the 
velocity  of  the  no-slab  case  should  be  recovered.  This  seems  to  be  the  case. 
However,  differences  in  the  temperature  profiles  persist  for  long  times. 

In  order  to  study  ablation  layer  Rayleigh-Taylor  instabilities,  numerical 
experiments  were  performed  using  the  code  FAST20. 

As  mentioned  previously,  we  initialized  the  system  using  the  steady 
state,  semi -analytical  solutions  of  Orens  for  our  starting  temperature, 
density,  and  pressure  profiles.  A dispersion  relation  is  obtained  by 
initially  perturbing  the  density  profile  at  its  peak.  The  first  ten  modes 
are  excited  with  equal  amplitudes  and  random  phases  corresponding  to  a total 
initial  density  perturbation  of  eight  percent. 

Using  the  new  diagnostic  described  previously* the  amplitudes  of  the 
first  three  modes  were  determined.  An  example  of  a plot  of  these  amplitudes 
is  given  in  Figure  4.  Growth  rates  were  determined  from  the  exponentially 
, growing  regions  of  such  curves. 

Figure  5 shows  the  dispersion  relation  obtained  from  four  different 
runs  as  labelled.  The  first  three  runs  were  performed  with  a 60  cell  by  20 
I cell  mesh.  For  these  runs,  a ly  resolution  in  the  x-direction  was  used  for 

the  non-stretched  cells.  The  y-resolution  for  the  three  runs  were  1,2,  and 
4y  respectively.  For  the  last  run,  a 120  cell  by  20  cell  mesh  was  employed. 
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Hence,  &x  was  reduced  to  0.5u.  This  increased  resolution  produced  minimal 
changes  in  the  growth  rates. 

The  shaded  area  in  Figure  5 corresponds  to  the  classical  Rayleigh- 
Taylor  dispersion  relation, 

y * / kg;  k * 2tr/X.  (9) 

The  width  of  the  region  is  due  to  the  variability  of  g,  the  acceleration. 

Although  we  have  a constant  laser  flux  of  1013  watts/cm2,  the  mass  of  the 
system  decreases  because  of  outflow  and,  hence,  the  acceleration  increases. 

Another  important  aspect  of  our  research  on  tne  ablation  layer  Rayleigh- 
Taylor  instability  involves  a study  of  its  dynamic  stabilization.  Boris 
proposed  that  by  oscillating  the  laser  intensity  at  some  given  frequency,  the 
growth  of  some  of  the  Rayleigh-Taylor  modes  could  be  inhibited.  The  fre- 
quency is  picked  in  order  to  damp  the  most  destructive  modes. 

Our  most  recent  dynamic  stabilization  studies  using  FAST2D  are  shown  in 
Figure  6.  Note  that  these  calculations  were  performed  on  an  older  version 
of  the  code.  Hence,  the  temperature  of  the  boundary  is  specified  rather 
than  the  laser  flux,  and  the  diagnostic  for  determining  the  mode  amplitudes 
is  not  as  clear  as  the  example  shown  in  Figure  4. 

The  upper  row  of  figures  displays  the  time  evolution  of  the  density 
with  no  dynamic  stabilization.  Note  that  by  7nsec  the  shell  is  eaten  through. 

For  the  dynamically  stabilized  case,  shown  in  the  bottom  row,  the  shell 
maintains  its  integrity  past  the  9nsec  point.  The  temperature  was  oscillated 
with  a frequency  of  1.25  cycles/nsec  and  an  amplitude  of  125  eV.  It  is  displayed 
in  the  lower  left  hand  corner  along  with  the  target  velocity.  Note  that  for 
the  stabilized  case  (solid  line)  almost  twice  the  target  velocity  is  achieved 
prior  to  shell  break-through.  In  the  upper  right  hand  corner,  a comparison 
of  the  mode  amplitudes  with  and  without  dynamic  stabilization  is  presented. 

The  crosses  refer  to  the  first  mode  in  the  stabilized  case,  whereas,  the 
circles  and  squares  represent  the  first  and  second  modes  respectively  for 
the  non-stabilized  run.  Again,  it  is  clear  that  the  stabilization  has 
inhibited  the  Rayleigh-Taylor  growth  rates. 
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Another  potential  cause  of  instabilities  on  the  surface  of  the  target 
shell  arises  from  the  use  of  thin  membranes  to  support  inner  shells  in 
multi-shell  fusion  targets.  These  support  membranes  (webs)  can  cause 
inhomogeneities  in  the  density  distribution.  We  have  studied  the  instabilities 
arising  from  these  in  a schematic  way. 

FAST2D  was  used  in  a slab  geometry  mode.  Besides  our  initial  density 
profile,  the  membrane  was  represented  by  an  osculating,  Gaussian  density 
perturbation  consisting  of  the  same  material  as  the  target. 

A uniform  amount  of  mass  for  each  y-row  was  distributed  between  two 
x-cells.  The  center  of  mass  of  the  membrane  lied  along  a prescribed  Gaussian. 
Because  of  grid  stretching  for  large  x,  and  finite  resolution  effects,  this 
uniform  membrane  appears  as  the  irregular  shape  in  the  upper  left  hand  corner 
of  Figure  7.  The  density  was  chosen  to  be  0.125  g/cc  in  the  unstretched 
cells.  For  comparison,  the  maximum  slab  density  is  0.315  g/cc.  The  rest  of 
Figure  7 gives  the  time  evolution  of  the  density.  One  sees  that  the  presence 
of  the  membrane  gives  rise  to  a large  "T-shaped"  perturbation  in  the  ablation 
region.  A detailed  study  of  the  early  timesteps  shows  that  this  mass  flow 
arises  from  thermal  conduction  effects.  Initially,  the  membrane  is  colder 
and  denser  than  the  surrounding  plasma.  The  membrane  heats  and  loses  mass, 
but  the  heating  is  much  more  rapid.  This  leads  to  a pressure  profile  that 
forces  mass  towards  the  center  in  the  y-direction  and  out  away  from  the  slab 
in  the  x-di recti  on.  Hence,  the  "T-shape"  is  formed. 

Further  tests  were  performed  using  a less  dense  membrane  and  one  dis- 
placed away  from  the  target.  The  upper  sequence  of  Figure  8 shows  how  the 
disturbance  develops  when  the  membrane  density  was  reduced  by  a factor  of 
five.  The  lower  half  is  for  the  case  when  the  membrane  is  displaced  5u 
away  from  the  target.  For  both  cases,  the  characteristic  "T-shape"  still 
develops,  albeit,  more  slowly. 

Figure  9 contains  the  first  three  mode  amplitudes  for  the  first  membrane 
test  discussed.  For  comparison,  the  first  mode  amplitude  from  our  standard 
density  perturbation  run  is  represented  by  the  straight  line.  One  sees  that 
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for  the  two  cases  the  growth  rates  are  comparable.  These  calculations  have 
given  a strong  Indication  that  these  support  webs  can  lead  to  large 
instabilities  in  the  ablation  process.  Further  work,  both  numerical  and 
experimental,  should  be  undertaken  to  study  this  in  greater  detail. 

All  high  Intensity  lasers  produce  non-uniform  intensity  patterns.  These 
intensity  variations  will  cause  non-uniform  heating  leading  to  inhomogenities 
in  the  ablation  process.  The  question  of  whether  these  variations  in  laser 
energy  are  transmitted  to  the  inner  edge  of  the  target  shell  is  an  important 
one. 


In  order  to  try  to  answer  this  question,  we  performed  some  calculations 
in  which  a 5%  sinusoidal  perturbation  of  the  laser  Intensity  was  introduced. 

In  Figure  10,  the  case  for  a 200y  wavelength  disturbance  is  presented.  As 
can  be  seen,  the  ablation  rate  variation  follows  the  laser  Intensity  varia- 
tion and  the  target  buckles.  At  least  for  this  case,  the  transverse  heat 
conduction  is  not  large  enough  to  smooth  out  the  thermal  inhomogeneities. 

We  also  undertook  a study  of  the  effects  of  varying  the  wavelength  of 
the  perturbation.  It  was  found  that  much  shorter  wavelength  disturbances 
did  not  appreciably  affect  the  target.  This  can  be  understood  from  the 
following  argument.  The  sinusoidal  disturbance  in  the  y-direction  should 
lead  to  a damped  thermal  disturbance  traveling  in  the  x-direction.  This 
disturbance  should  e-fold  over  a distance  comparable  to  the  wavelength  of  the 
perturbation.  Hence,  disturbances  of  wavelength  10y  will  not  appreciably 
affect  the  target  which  is  on  the  order  of  10-40y  away  from  the  absorption 
region,  whereas  the  200y  disturbance  strongly  warps  the  target.  These 
results  should  be  of  use  in  the  design  of  the  high  powered  lasers  necessary 
in  fusion  research. 

Two  different  parameterization  schemes  of  the  planetary  boundary  layer 
(PBL)  have  been  tested  and  proposed  for  incorporation  into  the  NRL  3D 
tropical  cyclone  model.  Table  1 shows  all  the  different  PBL  parameterization 
schemes  that  have  been  tested. 


PBL  PARAMETERIZATION  TESTED 


CONFIGURATION 

BASIC  MODELS 

VARIATIONS 

Single-Layer 

Aerodynamic  Drag 

Coefficient 

Constant  CD 

Wind  Dependent  Cq 

i 

Deardorff  (1972) 

Constant  PBL  Height 

ih  1,0 

Matching 

Arya's  Universal 

Functions 

Yamada's  Universal 

Functions 

Multi-Layer 

Constant  K Coefficient 

i 

1 

j 

Unconditional  Mixing 

Conditional  Mixing 

O'Brien's  Parabolic  K j 

Unconditional  Mixing 

- | 1— — . ..  . 

Conditional  Mixing 

I 

Busch  et  al^  Model 

I 

Various  Mixing  Lengths 

Meller  and  Yamada's  Level 

II  Closure  , 

i 

Two- Parameter  Model  i 

(Marchuk,  1977) 

I 

Ri  Adjustment  (Chang  1979) 

I 


The  first  type  is  to  be  used  when  there  will  be  several  layers  In  the 
tropical  cyclone  model.  For  this  purpose,  an  efficient,  multilayer  model 
for  resolving  the  thermal,  moisture,  and  momentum  fields  In  the  PBL  has  been 
developed.  The  surface  heat  flux  can  be  generated  by  two  mechanisms: 

(1)  the  convective  mixing  depending  on  the  temperature  difference  between 
the  surface  and  the  screen  level;  and  (2)  the  mechanical  mixing  depending  on 
the  wind  stress.  In  the  PBL,  the  heat  and  momentum  exchanges  are  computed 
by  a Richardson  number  (Ri)  adjustment  scheme.  Heat  and  momentum  exchanges  are 
computed  mainly  due  to  thermal  Instability  arisinq  under  convectively  unstable  con- 
tlons.  The  computational  scheme  is  based  on  the  concepts  of  conservation  of 
enthalpy  and  momentum.  Case  studies  show  good  agreement  between  model  results 
and  the  Great  Plains  Experiment  dataJ  Figure  11  shows  the  evolution  of  the 
PBL  as  illustrated  by  the  potential  temperature  profile.  The  model  is  capable 
of  reproducing  the  PBL  heights,  the  overshooting,  the  sharpening  of  the  inver- 
sion, the  slightly  stable  lapse  rate  in  upper  PBL,  and  the  unstable  lapse  rate 
in  lower  PBL  that  are  evident  in  the  observation. 

Since  no  detailed  PBL  observations  have  ever  been  attempted  during  tropical 
cyclones  over  the  ocean,  comparison  has  been  made  between  this  efficient  model 
and  the  more  sophisticated  and  much  more  expensive  PBL  model  of  Busch  et  al. 

The  two  models  are  Integrated  from  the  same  initial  condition  typical  of  a 
hurricane  PBL.  The  results  shown  in  Figure  12  after  12  hours  of  integration 
of  the  two  models  agree  remarkably  well  in  PBL  depth,  the  strength  and  height 
of  the  inversion,  the  radial  velocity,  the  mean  potential  temperature,  and  the 
mean  water  vapor  content.  The  computational  cost  of  the  Busch  et  al  model 
is  12  times  that  of  the  efficient  model. 

Since  no  explicit  diffusion  coefficient  is  involved  in  this  efficient 
model,  the  time  step  is  not  restricted  by  the  computational  stability 
requirement  of  the  diffusion  term.  The  five  minute  time  step  used  in 
obtaining  these  results  is  what  is  currently  allowed  by  the  split-explicit 
Integration  scheme  used  in  the  NRL  tropical  cyclone  model. 


Reported  in  Chang,  S.  W.,  1979:  "An  Efficient  Parameterization  of  Convective 
and  Non-Convectlve  PBL  for  Use  in  a Dynamic  Model",  to  appear  in  Journal 


The  second  boundary  layer  parameterization  is  designed  for  use  in  the 
NRt  tropical  cyclone  model  in  which  the  whole  boundary  layer  Is  represented 
by  one  model  layer.  Here,  the  generalized  similarity  theory  is  utilized. 

The  "matching  technique"  discussed  by  Blackadar  and  Tennekes  is  used  to 
obtain  a set  of  PBl  equations  that  has  the  characteristics  of  both  the  mixed 
layer  and  the  constant  flux  layer.  By  scaling  this  matched  set  of  PBL 
equations,  non-dimensional  gradient  equations  containing  universal  functions 
A,  B,  and  C result.  By  using  the  proper  scales  and  observational  data, 

Yamada  computed  the  values  of  these  universal  functions  that  show  the  least 
scattering.  The  polynomial  form  of  Yamada's  universal  functions  also  facili- 
tates vectorization  of  the  computer  code. 

To  further  test  Yamada's  functions  under  extreme  conditions,  the  matched 
equations  are  used  to  diagnose  the  surface  fluxes  based  on  the  interior  flow 
characteristics  of  twenty-four  different  12  hour  integrations  of  the  Busch 
et  al  model.  The  conditions  under  which  the  Busch  et  al  model  are  integrated 
are  highly  non-steady  with  stabilities  varying  from  the  strongly  stable  to 
the  strongly  unstable.  The  diagnosis  using  Yamada's  functions  yields  a high 
bias  in  both  the  unstable  and  stable  cases  when  evaluating  the  surface  stress 
as  shown  by  Figure  13.  The  diagnosis  produces  a high  bias  in  the  unstable 
cases  and  a low  bias  in  the  stable  cases  in  evaluating  the  surface  water  vapor 
flux.  Very  good  agreement  is  found  in  surface  water  vapor  flux  for  nearly 
neutral  cases  (Figure  14).  These  biases,  however,  decrease  with  time  when  the 
PBL  is  approaching  quasi-steady  state. 

The  PBL  parameterization  using  the  matching  technique  has  been  tested  with 
a 20  and  3D  tropical  cyclone  model.  Preliminary  results  show  that  cyclones 
with  very  good  dynamic  and  thermodynamic  structure  are  generated  and  maintained 
at  realistic  intensity. 
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IV.  CONCLUSION 


In  sections  II  and  III,  we  described  many  of  the  specific  computational 
advancements  and  scientific  studies  that  have  been  completed  under  this 
contract.  Since  many  of  the  codes  mentioned,  and  especially  subsections  of 
them,  are  used  in  other  NRL  projects,  these  accomplishments  have  other 
uses  that  we  have  not  described.  For  example,  some  of  the  modifications  in 
the  ID  code  have  been  incorporated  into  the  code  FLAME  which  is  used  in 
reactive  flow  studies  involving  complicated  chemical  energy  release.  Also, 
the  2D  tests  involving  varying  grid  resolution  and  the  improved  diagnostic 
for  determining  ablation  layer  irregularities  are  of  use  in  the  calibration 
studies  of  turbulence  closure  models. 
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FIGURE  CAPTIONS 

1.  A comparison  of  the  numerical  results  of  DIFFU1  (points)  with  the 
analytical  solutions  (lines)  for  the  thermal  diffusion  test.  On  the  left, 
the  "shock"  position  vs.  time  is  shown  and  on  the  right,  the  temperature 
profiles  are  displayed. 

2.  The  thermal  burn  through  of  1 p slabs  of  different  materials  as  labeled  is 
modeled  using  DIFFU1.  The  decrease  of  the  thermal  conductivity  through 
the  slabs  is  evident. 

3.  A comparison  of  thin  shell  burn  through  for  1 and  2 p slabs  and  a no-slab 
control  run  for  longer  times.  Although  the  velocities  converge  at  longer  dis- 
tance from  the  slab,  the  differences  in  temperature  and  "shock"  position 
persist. 

4.  The  first  three  mode  amplitudes  for  one  of  the  Rayleigh-Taylor  simulation 
runs  are  displayed.  The  growth  rates  can  be  determined  accurately  in  the 
exponentially  rising  region  (here  ~ 3-5  nsec). 

* I -I 

i -■ 

5.  A comparison  of  all  of  the  determined  growth  rates  with  the  classical 
Rayleigh-Taylor  dispersion  relation  (shaded  area).  Note  that  the  improved 
resolution  does  not  significantly  change  the  growth  rates. 

6.  A demonstration  of  the  dynamic  stabilization  of  some  of  the  most  damaging 
Rayleigh-Taylor  unstable  modes.  The  upper  sequence  shows  a time  evolution 
of  the  slab  when  an  initial  density  perturbation  is  applied  at  the  density 
peak.  The  lower  sequence  is  for  the  same  case  but  with  dynamic  stabiliza- 
tion. As  can  be  seen,  the  shell  maintains  its  integrity  longer  when 
stabilized.  The  lower  left  hand  corner  contains  a plot  of  the  driving 
temperature  for  the  stabilized  case  (upper  half)  and  a comparison  of 
target  velocity  for  the  stabilized  (solid  line)  and  non-stabilized  (dashed 
line)  calculations  (lower  half).  Note  that  almost  twice  as  large  a target 
velocity  Is  achieved  in  the  stabilized  case.  In  the  upper  right  hand 
corner  Is  a comparison  of  the  growth  rates  for  the  two  cases.  The 
crosses  represent  the  first  mode  amplitude  for  the  stabilized  case  and 
the  circles  and  squares  show  the  first  and  second  modes,  respectively,  for 
the  non-stabilized  calculation. 
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7.  The  time  evolution  of  the  ablation  of  a target  shell  in  the  presence 
of  a support  web. 

8.  The  same  as  Figure  7 except  that  in  the  upper  sequence  the  membrane's 
density  has  been  reduced  by  a factor  of  five,  and  in  the  lower  one  the 
membrane  has  been  displaced  5/i  away  from  the  target. 

9.  The  first  three  mode  amplitudes  for  the  web  calculation  that  was  shown 
in  Figure  7.  For  comparison,  the  first  mode  amplitude  of  the  standard 
Rayleigh-Taylor  test  is  represented  by  the  solid  line. 

10.  The  effect  of  a 5%  variation  in  the  laser  intensity  on  a plasma  target 
Is  shown. 

11.  Potential  temperature  profiles  at  different  local  mean  times  computed  by 
the  Richardson  number  adjustment  model . 

12.  Comparison  of  the  Richardson  number  adjustment  model  with  the  Busch  et  al 
model.  Results  shown  are  after  12  hours  of  simulation  time.  The  two 
models  are  initialized  with  the  same  initial  conditions  typical  of  a 
tropical  cyclone. 

13.  Comparison  of  surface  stresses  as  computed  by  the  Busch  et  al_  model  and 
diagnosed  by  the  generalized  similarity  theory  with  Yamada’s  functions. 

14.  Same  as  figure  13  except  for  surface  water  vapor  flux. 
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ASYMPTOTIC  RESULTS  FOR 
fHIN  SHELL  BURN  THROUGH 


A Ly  = 20fj 


ABLATION  LAYER  RAYLEIGH  - TAYLOR  INSTABILITY 


RAYLEIGH-TAYLOR  INSTABILITIES 
GENERATED  BY  A SUPPORT  WEB 


Figure  9 


5%  VARIATION  IN  LASER  INTENSITY 
FOR  1013  w/cm2  IRRADIATED  D-T 


Figure  10 
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